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After Anderson’s prediction of disorder-induced insulating behavior, extensive work found no 
singularities in the density of states of localized systems. However, Johri and BliatlP! recently 
uncovered the existence of a non-analyticity in the density of states near the band edge of non¬ 
interacting systems with bounded disorder, in an energy range outside that captured by previous 
work. Moreover, this feature marks the boundary of an energy range in which the localization is 
sharply suppressed. Given strong current interest in the effect of interactions on disordered systems, 
we explore here the effect of a Hubbard U interaction on this behavior. We find that in ensembles 
of small systems a cusp in the density of states persists and continues to be associated with a sharp 
suppression of the localization. We explore the origins of this behavior and discuss its connection 
with many-body localization. 


I. INTRODUCTION 

Unconventional superconductivity and other novel 
electronic behaviors in transition metal oxides continue 
to drive the study of strong electron correlations. Dop¬ 
ing which breaks translational invariance is generally re¬ 
quired to obtain these novel behaviors, raising the ques¬ 
tion of how disorder affects strongly correlated systems. 
Meanwhile, since the seminal work of Anderson in 1958,® 
localization of single-particle states by disorder has been 
the subject of intense study in a wide range of systems.® 
Experiments have suggested that electron-electron inter¬ 
actions might drive a delocalization transition,® and there 
has been recent rapid progress in addressing the question 
of how interactions affect localized systems.® 1 ® 

In this paper we ask how electron-electron interac¬ 
tions affect several novel features re cen tly pointed out 
in non-interacting disordered systems? 2 After Anderson 
first predicted localization, significant effort went toward 
establishing whether or not any feature in the density 
of states was associated with localization. Edwards and 
Thoulesd 1 ^ showed that the density of states is analytic 
in a wide energy range around the middle of the band for 
all continuous distributions, and Wegner 12 showed that 
the density of states neither vanishes nor diverges any¬ 
where inside the band for Gaussian distributed disorder. 
It therefore came as a surprise when Johri and BhatlP ® 
uncovered the existence of a non-analyticity in the den¬ 
sity of states of systems with bounded disorder. Their 
result does not contradict earlier work. Rather, the non- 
analyticity they found occurs outside the energy range 
addressed by Edwards and Thouless. 

Beyond the novelty of this singularity, Johri and Blratt 
showed that a corresponding singularity exists in the in¬ 
verse participation ratio, a measure of localization, and 
demonstrated that these singularities mark the boundary 
of a region at the edge of the band dominated by resonant 
states, also known as Lifshitz states. Moreover, while this 
behavior persists in large systems and in two and three 
dimensions,® Johri and Bhatt demonstrated that the key 


features can be seen and understood analytically in the 
simple case of an ensemble of two-site systems.® 

Here we consider an ensemble of two-site Anderson- 
Hubbard systems and calculate both the density of states 
and a generalized form of the inverse participation ratio 
applicable to interacting systems. Fig. [TJsummarizes our 
key results. The presence of local interactions results in 
a richer energy dependence of both quantities. In partic¬ 
ular, the minimum in the inverse participation ratio no 
longer occurs at the band edge but instead moves closer 
to the Fermi level. We present below an exploration of 
how these changes arise, highlighting two distinct types 
of resonance which can occur in interacting systems and 
their role in our results. 


II. THE MODEL 

The Anderson-Hubbard model is a tight-binding 
Hamiltonian which includes nearest-neighbor hopping t , 
on-site Coulomb repulsion U, and site potentials chosen 
from a uniform distribution of width W. 

7-L = t ^2 cl^jy + ^ Uh^hii + y^(ej - n)h ia (l) 

i i,<r 

c\ a is the creation operator for an electron with spin a 
at lattice site i, fii a = c\ a Ci a , and (■ i,j) refers to nearest 
neighbor pairs. The site potentials €i are chosen from 
the distribution P(ei) = 0(W/2 — ]e»|)/W where 0 is 
the Heaviside function, y, is the chemical potential. 

We consider here ensembles of two-site systems so 
i,j = 1,2, and we focus on the case of half filling so 
H = U /2. For each system in an ensemble, the site po¬ 
tentials ei and ti are different but all are chosen from 
the same distribution. Following Johri and Bhatt, it is 
convenient to use the coordinates 

x = 61 ~t- C2 and y = £l 62 . (2) 
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FIG. 1: (Color online) Panels (a)-(d) show the ensemble-averaged density of states for W/t = 12 and U/t = 0, 2, 4 and 
8 respectively. The thin black lines indicate the atomic-limit density of states. Panels (e)-(h) show the ensemble-averaged 
generalized inverse participation ratio at the same U values. The vertical dashed lines mark energies at which a cusp appears 
in the density of states and a sharp drop occurs in the generalized inverse participation ratio. Ensembles consist of 100-150 
million systems and the energy bin width is 0.019t-0.025t. 
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FIG. 2: (Color online) In the main panel, each point in the 
diamond specifies the site potentials for one 2-site system in 
an ensemble with disorder strength W/t = 12. Regions are 
labeled according to the number of particles in the ground 
state for interaction strength U/t = 8. Insets show typical 
single-system density of states plots for each region with the 
transition corresponding to each peak labeled. Energy bin 
width 0.02t 


The phase space of all systems in an ensemble is repre¬ 
sented by the diamond in the main panel of Fig. [2j 


III. THE DENSITY OF STATES 


The density of single-particle states for a single system 
is the average of the local density of states on the two 
sites calculated from the Lehmann representation of the 


retarded Green’s function (see appendix). 

pH = * (piM + p 2 (w)) (3) 

where Pi(u>) = — — Im G^(uj) (4) 

7r 

The result consists of a small number of ^-functions cor¬ 
responding to the allowed single-particle transitions from 
the ground state. 

The main panel of Fig. [2] shows the number of par¬ 
ticles in the ground state in each region of the phase 
space, and the insets show typical density of states plots 
in each of these regions. For systems with a 0-particle 
ground state, there are two distinct peaks corresponding 
to transitions to the 1-particle bonding (lb) excited state 
or the 1-particle anti-bonding (la) excited state. For 
convenience we label the transition from the ?r-particle 
ground state to excited state i the n:i transition. In this 
notation these two peaks correspond to the 0:1b and 0:1a 
transitions. Each includes a spin-up and a spin-down 
contribution which, because our Hamiltonian is param¬ 
agnetic, occur at the same energy. For systems with a 
1-particle ground state, there are five distinct peaks cor¬ 
responding to transitions to the 0-particle state and to six 
possible 2-particle states. Three of the 2-particle states, 
the triplet states (2t), are degenerate. The remaining 2- 
particle states (2a, 2b, and 2c) are linear combinations 
of the singlet state |s) = (| t4-) — I Y())/V^ an d the 
double-occupancy states |20) and |02) in the Fock ba¬ 
sis. For systems with a 2-particle ground state, there are 
four distinct peaks corresponding to transitions to two 
1-particle states (la and lb) and two 3-particle states 
(3a and 3b). The transition options from the 3-particle 
ground state mirror those of the 1-particle ground state, 
and similarly the 4-particle options mirror those of the 
0-particle ground state. 

The density of states of an ensemble of two-site systems 
simply averages the contributions of all the systems. The 
density of states plots in Fig. [l] (a)-(d) are histograms 
of the contributions in a sequence of frequency bins re- 
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suiting in a smooth density of states for a large but finite 
number of systems. 

In considering the evolution of the DOS with interac¬ 
tion strength, the atomic limit (t = 0) provides a useful 
point of comparison. This is indicated by the narrow 
black lines in panels (a)-(d) of Fig. |T] Without interac¬ 
tions all sites are either empty or doubly occupied, and 
the atomic-limit DOS is identical to the distribution of 
site potentials: a, plateau of height 1/IF. On-site interac¬ 
tion results in singly-occupied sites which at U = 2 create 
a narrow raised region at the band centre. The width of 
the contribution from singly occupied sites grows with U, 
matching that from empty and doubly occupied sites at 
U = 4 and exceeding it at U = 8. 

With hopping, the occupancy of each site is no longer 
integer, and there is an increase both in the number of 
distinct energies at which transitions can occur and in 
the range of weights associated with these transitions. In 
the non-interacting case, this results in the two peaks at 
the band edge highlighted by Johri and BliattP-® With 
interactions, the DOS gains even more structure. The 
zero-bias anomaly in panel (d) has been explored in detail 

elsewhere! 13 ! 14 ! 

The particular feature on which we focus here is the 
persistence of sharp peaks, marked in Fig. |T] by verti¬ 
cal dashed lines. These peaks closely follow the edge of 
the contributions from empty and doubly occupied sites 
in the atomic limit, which might suggest that the cusps 
would be lost when U > W and all sites are singly oc¬ 
cupied. In fact, when U > W the shoulders seen around 
ui/t = ±9 in panel (d) evolve into sharp cusps, such that 
the DOS of the Mott insulator (not shown) looks like two 
copies of the U = 0 case above. 

To provide further insight into the origin of these 
peaks, Fig. [3] shows the DOS contributions from indi¬ 
vidual transitions. The first column indicates the transi¬ 
tion, and the second shows the DOS contribution. In the 
third column the phase space is shown with a blue curve 
which indicates all systems for which the specified transi¬ 
tion occurs at the energy marked with a vertical dashed 
line in column 2. As the energy E is increased, the shape 
of this curve for a given transition does not change but 
the curve shifts to the right. In the non-interacting case, 
the magnitude of the DOS contribution from each tran¬ 
sition in each system is the same. Therefore the overall 
magnitude of the DOS at a given energy is simply pro¬ 
portional to the length of the constant energy curve at 
that energy. In the interacting case, there is variation in 
the magnitudes of the DOS contributions from different 
transitions, meaning that the line length alone is not the 
sole determinant of the DOS magnitude at each energy. 
However, the main qualitative features of the DOS con¬ 
tribution from each transition are still reflected in the 
changes in line length. 

Consider, for example, the 1:2a transition. The con¬ 
stant energy curve has the same shape as the boundary 
between the 2-particle and the 1-particle ground state 
regions. At E = 0, the constant energy curve coincides 
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FIG. 3: (Color online) Transition specific information for an 
ensemble of two-site systems with W/t = 12 and U/t = 8. 
Each row corresponds to a different transition as labeled in 
column 1. Only transitions which contribute above the Fermi 
level are shown. The 1:2c transition is not shown because its 
contribution is negligible. Column 2 shows the contribution 
to the density of states made by this transition alone (red) 
as well as the full density of states (black) for comparison. 
Column 3 shows the constant energy curve (blue) in phase 
space at the energy marked with a dashed line in column 2. 
Column 4 shows the generalized inverse participation ratio 
for this transition at a set of points along the constant energy 
curve shown in column 3. 


with this boundary, and the transition makes its maxi¬ 
mum contribution to the DOS. As the energy is increased, 
the constant energy curve moves to the right declining in 
length roughly linearly with energy. In this way the ori¬ 
gin of the energy dependence of the DOS contributions of 
each transition can be seen, and from these the structure 
of the DOS as a whole. 

In particular, what do we learn from Fig. [3] regarding 
the origin of the DOS cusps? As in the non-interacting 
casef 2 the peaks appear where a constant energy curve 
sits nearly tangent to the edges of the phase-space di¬ 
amond, maximizing its length. (See the l:2t transition 
in Fig. HP Unlike the non-interacting case, many transi¬ 
tions contribute at the same energy, resulting in variation 








































4 


in the shape of the peaks. 


IV. THE GENERALIZED INVERSE 
PARTICIPATION RATIO 


To explore whether there is a connection between the 
structure in the density of states and variations in the 
strength of localization, we have also studied a gener¬ 
alization of the inverse participation ratio, a standard 
measure of localization. In non-interacting systems, the 
inverse participation ratio at an energy E is proportional 
to one over the size of a single-particle state with en¬ 
ergy E. In interacting systems, single-particle states are 
not well defined. Nonetheless, the following generalized 
inverse participation ratio (GIPR) can be defined: 


i M 


EiP?M 

E iPM? 


( 5 ) 


Pi(u>), the local density of states as defined above in 
Eqn. Q, has the form of a sum of weighted d-functions: 
Piioj) = ~ Et) where E t are the energies 

of single-particle transitions available from the ground 
state. In many numerical calculations only the local den¬ 
sity of states is known and only up to some nonzero en¬ 
ergy resolution. In such cases, when interactions are set 
to zero the GIPR has the same scaling behavior as the 
IPR in the limits of zero and infinite disorder, but the 
correspondence for intermediate disorder is less clear. 15 
In our case, however, we know the many-body eigen¬ 
states and can therefore distinguish the local-density-of- 
states contributions of individual transitions no matter 
how close they may be in energy. Regularizing the 5- 
functions in the LDOS would result in a loss of informa¬ 
tion. To avoid this we interpret Eqn. © as 

V w 2 - 

m= Jr E for u = Et (6) 

E* w ti\ 

and zero otherwise. In this way, in the absence of interac¬ 
tions the GIPR reduces to the IPR, expressed in terms of 
the single-particle eigenstate as JT |t/’?:| 4 / [Ei IV'il 2 ] • 
Note that in an interacting system, the GIPR is asso¬ 
ciated with the size of a transition between two many- 
body states. Column 4 of Fig. [3] shows the GIPR values 
for specific transitions in selected individual systems. For 
example, consider the 2:3a transition. The blue line in 
column 3 shows all systems for which the energy of this 
transition is 8.1 1. The system at y = 5 has a GIPR value 
of 1 for this transition. In this system, the LDOS associ¬ 
ated with this transition is large on one site and close to 
zero on the other. This is a transition which is localized 
primarily on a single site. Alternatively, the system at 
y = 0 has a GIPR of 0.5 for this transition. In this case, 
the LDOS associate with the transition is roughly equal 
on the two sites, and the transition is as extended as it 
can be in a two-site system. 


Column 4 of Fig. [3] shows GIPR values for each tran¬ 
sition for a selection of systems at different values of y 
along the constant energy curve shown in column 3. The 
GIPR values are independent of x. For most transitions 
the GIPR value increases as y increases. For example, 
2:3a transitions have GIPR values of 0.5 at y = 0 and 
these values increase rapidly to 1 at higher values of y. 

The 2:3b transition is particularly interesting: There 
are two minima in the GIPR, one at y = 0 and one at 
y = ±U/y/2. These two values represent the two distinct 
types of resonance which can occur in a strongly corre¬ 
lated system. The y = 0 resonance is the one present in 
non-interacting systems: The potentials on the two sites 
are the same. The y = U/y/2 resonance is unique to 
strongly correlated systems: Here the potential on one 
site is greater by U than the potential on the other site. 

From column 4 of Fig. [3] it is clear that the transitions 
contributing to the DOS at any given energy may have 
a wide range of GIPR values, both because of contribu¬ 
tions at the same energy from different transitions and 
because of variation in the GIPR within a single transi¬ 
tion. To obtain a single value at each energy we average 
over all the contributions. In the non-interacting case 
studied by Johri & BhatlP^, the magnitude of the DOS 
contribution from each transition is the same, so the IPR 
values can be directly averaged. In the interacting case 
we study here, there is wide variation in the magnitude of 
the DOS contributions. We therefore average the GIPR 
values weighted by the magnitude of the corresponding 
DOS contributions. Consider a transition t in system s 
at frequency u> s t- The DOS contribution for this tran¬ 
sition is p st and the GIPR value is G s t ■ We define the 
ensemble-average GIPR value 


(G M> 


Es Et G st (uj)pst(uj)5(uj - UJ st ) 

I2sI2tPst(u)6(u]-U) s t) 


Fig. a (e)-(h) show the ensemble-average GIPR versus 
energy for four values of U. At U = 0 there is a gentle 
upward curvature in the centre of the band mirroring the 
gentle downward curvature in the DOS, consistent with 
the usual picture that localization is stronger when the 
DOS is smaller. Then, after reaching a sharp maximum 
there is an abrupt drop in the IPR at the edges of the 
band. This drop has been showrP to come specifically 
from systems with extreme x values but y near zero. The 
states in these systems have been referred to as resonant 
or Lifshitz states. 

When U is nonzero, the GIPR, like the DOS, has much 
more structure. One pattern is an overall lowering of the 
GIPR over most of the width of the band as U is in¬ 
creased. This is consistent with weak interactions pro¬ 
viding screening which reduces localization. 

A second pattern is the persistence of region of sharply 
suppressed localization. In particular this region moves 
towards the band center as U increases. What causes 
this abrupt drop? Changes in the GIPR can occur at fre¬ 
quencies at which a particular transition starts or stops 
making DOS contributions or due to variations in the 
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GIPR value within a specific transition. The latter effect 
is generally smoother and smaller in magnitude. The 
sharp drops in GIPR which move toward the center of 
the band with increasing U are due to equal-site-potential 
resonances, as in the noninteracting case. There do not 
appear to be such dramatic features associated specifi¬ 
cally with the y = U/y/2 resonances unique to strongly 
correlated systems. 

V. DISCUSSION 


ing mostly on spin systems. It can be shown that the 
Hamiltonian of any system with a Hilbert spac e of di¬ 
mension 2 n can be expressed in terms of spinors E2I22I In 
this spinor picture, eigenstates differ in the flipping of 
individual spinors. Therefore, the GIPR presented here, 
which is a measure of the localization of transitions be¬ 
tween many-body eigenstates, is providing information 
on the localization of these spinors. 

Moving away from half filling and exploring larger sys¬ 
tems will be important directions for further study. 
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Appendix: The Green’s function 


Gu( u ) — ~ (Guff (w) + G^j_(w)) (A.l) 


where 


GLa M = E 


\(j’n\cl a \^o )\ 2 
uj — (fl n — H 0 ) + iy 

\(ipn\c ia \i>o )\ 2 

w + - H 0 ) + iy 


(A.2) 


Because the Hamiltonian conserves spin, G^ a/3 (uj) is zero 
for a 7^ f3. ipo is the many-body eigenstate with the low¬ 
est grand potential fl 0 = E 0 — y,N 0 where E 0 and N 0 
are the energy and the particle number corresponding to 
this state. i /> n are all other many-body eigenstates with 
grand potentials Q n . In non-interacting systems the dis¬ 
tribution of this local Green’s funct ion is an order pa¬ 
rameter for the Anderson transition ) 21 ! 22 ! and it is also 
being explored as a measure of the many-body localiza¬ 
tion transition!^! 
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